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ABSTRACT 

We present two-dimensional high-resolution hydrodynamical simulations in spher- 
ical polar coordinates of a Type la supernova interacting with a constant density in- 
terstellar medium. The ejecta are assumed to be freely expanding with an exponential 
density profile. The interaction gives rise to a double-shocked structure susceptible 
to hydrodynamic instabilities. The Rayleigh- Taylor instability initially grows, but the 
Kelvin-Helmholtz instability takes over, producing vortex rings. Provided the simula- 
tion is initiated early in the evolution with a perturbation > 1%, the instability reaches 
its fully developed nonlinear strength within 5 doubling times. The further nonlinear 
evolution does not depend on the initial conditions. Considering the small initial radii 
of Type la supernovae, they are likely to reach this fully developed phase. The nonlinear 
instability initially evolves toward longer wavelengths and eventually fades away when 
the reverse shock front is in the flatter part of the supernova density distribution. Based 
on observations of X-ray knots and the protrusion in the southeast outline of Tycho's 
supernova remnant, we include clumping in the ejecta because these features cannot 
be explained by instabilities growing from small perturbations. The clump interaction 
with the reverse shock induces Rayleigh- Taylor and Kelvin-Helmholtz instabilities on 
the clump surface that facilitate fragmentation. In order to survive crushing and to 
have a bulging effect on the forward shock, the clump's initial density ratio to the sur- 
rounding ejecta must be at least 100 for the conditions in Tycho's remnant. The 56 Ni 
bubble effect may be important for the development of dumpiness in the ejecta. The 
observed presence of an Fe clump would then require a non-radioactive origin for this 
Fe, possibly 54 Fe. The large radial distance of the X-ray emitting Si and S ejecta from 
the remnant center indicates that they were initially in clumps. 

Subject headings: hydrodynamics — instability — supernova remnants — supernovae: 
general — supernovae: individual (SN 1572) 
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1. INTRODUCTION 

Type la supernovae (SNe) are thought to arise from the thermonuclear explosion of a white 
dwarf that accretes mass from its binary companion (Whelan Sz Iben 1973; Chevalier 1981a; 
Nomoto, Thielemann, Sz Yokoi 1984; Harkness 1991). Both Chandrasekhar and sub-Chandrasekhar 
mass models have been considered for the mass of the white dwarf at the time of the explosion. 
In the sub-Chandrasekhar mass models, He detonation occurs on the surface of the white dwarf 
(Woosley & Weaver 1994; Livne & Arnett 1995). Although these models have some attractive fea- 
tures, they seem less able than Chandrasekhar mass models to explain the light curves and spectra 
of Type la supernovae (e.g., Hoflich & Khokhlov 1996). In Chandrasekhar mass models, ignition is 
initiated in the high density central regions of the white dwarf. In one type of model, the burning 
front propagates subsonically through the white dwarf. A popular example of this deflagration 
model is the W7 model of Nomoto et al. (1984), which can approximately reproduce the observed 
light curves. Branch et al. (1985) found that the model can reproduce the observed spectra if there 
is mixing of matter with a velocity > 8,000 km s _1 , although Harkness (1991), in more detailed 
models, found that mixing may not be necessary. In another class of Chandrasekhar mass models, 
the initial subsonic burning front, or deflagration, makes a transition to a supersonic burning front, 
or detonation. This is the delayed detonation, or DD, model (Khokhlov 1991). 

Although observations of extragalactic Type la supernovae have given insights into the nature 
of the explosions, a clear determination of the supernova model has not been possible. Observa- 
tions of Galactic supernova remnants provide another window on the structure and composition of 
supernovae. Likely remnants of Type la explosions are Tycho's supernova and SN 1006. Tycho's su- 
pernova (SN 1572) appears to be a Type la supernova based on its reconstructed light curve (Baade 
1945), although it may have been subluminous (van den Bergh 1993). Schaefer (1996) noted that 
the light curve allows a Type lb origin, but the lack of evidence for stellar wind or photoionization 
effects in the surroundings of Tycho's supernova argues against a massive star origin. At X-ray and 
radio wavelengths, Tycho shows a similar morphology, characterized by a limb-brightened circular 
shell. The outer edge is sharply defined and is interpreted as the forward shock front. Seward, 
Gorenstein, & Tucker (1983) found evidence for 400 X-ray clumps in the Einstein observatory im- 
age. At optical wavelengths, an Ha emission filament, which is thought to be from a fast shock 
moving into a surrounding neutral material, is observed along the NE perimeter, supporting the 
shock interpretation. The remnant outline spans 8' in diameter, corresponding to a linear radius 
of ~ 3 pc for an assumed distance of 2.5 kpc (Raymond 1984). The HI absorption show that the 
environment surrounding Tycho is fairly uniform, without a prominent density gradient (Reynoso 
et al. 1999). This property makes it especially suitable for detailed study. 

Hamilton et al. (1986) calculated a one-dimensional hydro dynamical model to study the X- 
ray spectrum of Tycho's remnant and obtained good agreement with spectra from the Einstein 
Observatory and other X-ray observatories. They took a three-layer composition structure and 
were successful in reproducing the observed Si line emission. However, they chose a constant 
density profile for the supernova. When this type of supernova density profile interacts with a 
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constant density ambient medium, the outer shocked layers (where Si was placed in the model) 
end up dense and cool. A constant density profile is probably a poor representation of an exploded 
Type la supernova. Dwarkadas <fe Chevalier (1998, hereafter DC98) recently reviewed the types 
of explosion models that reproduce the basic properties of Type la supernovae and found that an 
exponential profile is the best approximation for the density profile overall, although there can 
be significant deviations from this profile. When such a profile interacts with a constant density 
medium, the shocked supernova gas has an approximately constant temperature. In order to 
reproduce the observed spectral properties of Tycho's remnant, DC98 suggested that either the 
supernova gas has dense clumps or there was an early phase of interaction with a circumstellar 
medium. This was on the basis of one-dimensional hydrodynamic models for the interaction of a 
Type la supernova with a surrounding medium. 

In order to make further progress, we extend the hydrodynamic modeling effort to two dimen- 
sions. One motivation for this is the investigation of instabilities that result from the deceleration of 
the supernova ejecta by the ambient medium. Chevalier, Blondin, & Emmering (1992) previously 
investigated the interaction of a power law supernova density profile with an ambient medium. In 
this case, the flow approaches a self-similar state, although turbulent motions are continuously fed 
by the instability. In the case of an exponential profile, the flow is no longer self-similar; the density 
profile effectively evolves from a steep power law profile to a flatter one. Chevalier et al. (1992) 
investigated the instability for a range of power law indices and found that the flow remained qual- 
itatively similar. In view of this, we expect the exponential case to be qualitatively similar while 
the steep part of the density profile is interacting with the surroundings, even though the flow is 
no longer self-similar. Our aim is to determine the nature of the evolution. Dwarkadas (2000) has 
already made some numerical studies of this phenomenon. 

Another part of our effort is to include density inhomogeneities in the freely expanding super- 
nova gas. The ASCA observatory has made possible the imaging of Tycho's remnant in moderately 
narrow X-ray bands (Hwang & Gotthelf 1997). There are two X-ray emitting knots on the eastern 
side of the remnant where there is a protrusion in the remnant's outline (Vancura et al. 1995; 
Hughes 1997; Hwang & Gotthelf 1997). One of these knots is strong in Si line emission and the 
other is strong in Fe line emission. With the ROSAT HRI, Hughes (1997) has measured the proper 
motion of features in the X-ray emission. He found that these two knots show approximately unde- 
celerated motion. On the other hand, radio observations over the whole remnant show substantial 
deceleration, with Vt/R~ 0.47 (Reynoso et al. 1997 and references therein), where V is the shock 
velocity and R is the radius. The implication is that the knots represent not only an inhomogeneous 
composition, but also higher density clumps that have not been decelerated. At a distance of 2.5 
kpc, the velocities of the knots are about 8,300 km s~ x (Hughes 1997). 

The other well observed Type la supernova remnant, SN 1006, does not show evidence for the 
expansion of X-ray emitting clumps. However, absorption line studies of the Schweizer-Middleditch 
star behind the remnant have provided a powerful probe of the supernova remnant structure and 
composition (Wu et al. 1983, 1997; Fesen et al. 1988). These observations show evidence for central 
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cool Fe II that is compatible with expectations for the unshocked Fe rich matter in the central region 
of a Type la supernova, although the amount of Fe is smaller than expected (Hamilton et al. 1997 
and references therein). In addition, there is evidence for Si II absorption that is redshifted over 
the velocity range 3, 000 - 7, 000 km s" 1 (Wu et al. 1997; Hamilton et al 1997). Fesen et al. (1988) 
initially argued that the Si rich matter is in a clump on the far side of SN 1006. However, Hamilton 
et al. (1997) more recently argued for a spherically symmetric model for the supernova ejecta in 
which there is a transition from Fe to Si rich gas at a velocity of 5, 600 km s" 1 and the supernova 
remnant is considerably more extended on the back side because of a low interstellar density in that 
direction. They argue against a Si clump, but appear to assume that the clump gas would have 
been shocked and cooled. Another possibility is that low ionization absorption is due to unshocked 
matter in the clump. The clump may be moving at about 7, 000 km s _1 and matter stripped from 
the clump may lead to the matter extending to lower velocities. The remnant can then have a 
shape that is closer to spherical symmetry; the presence of a clump implies that the ejecta are not 
in spherical shells, but this is also indicated in Tycho's remnant. 

Our aim here is to examine how inhomogeneities in the ejecta interact with the decelerated 
interaction region. The observed properties can place constraints on the nature of the inhomo- 
geneities, which can be important for the explosion models. 

The plan of our paper is as follows. In § 2, we discuss the density structure of a Type la 
supernova and consider mechanisms that might give rise to inhomogeneities in the ejecta. Our 
method is given in § 3. Results on the growth of instabilities from small perturbations are in 
§ 4. The evolution of nonlinear inhomogeneities is given in § 5. Our results are discussed in the 
context of observations of Tycho's supernova remnant in § 6. Implications for models for Type la 
super novae are in § 7 and the conclusions are in § 8. 



2. DENSITY STRUCTURE IN TYPE la SUPERNOVAE 

Soon after a supernova explosion, the ejecta are expected to be freely expanding, with velocity 
v = r/t. In this phase, DC98 found that an exponential density profile generally describes the 
density distribution obtained from numerical 1-D (one-dimensional) explosion models, including 
delayed detonation, pulsating delayed detonation, He detonation, and deflagration models. The 
exponential density profile is given by 

Psn = Aexp(-v/v e ) t~ 3 , (1) 

where A is a constant and v e is the velocity scale height, which is determined by the total ejecta 
mass M and explosion energy E. This expression allows for spherical expansion; an element of 
supernova gas moves with a constant velocity and has its density drop as t -3 . Integrating the 
density and energy density over space and time gives (DC98) 



M = AttAvI, E = A&nAv b e , 



(2) 
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or 
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where M c ^ = 1.4M is the Chandrasekhar mass and E51 is the explosion energy in units of 10 51 
ergs. A larger ratio of E/M gives a larger velocity scale height and produces a flatter density profile 
at a given velocity. 

In some Type la models and core collapse models, the outer density distribution has been 
described by a power law. Then, the ejecta density is given by p oc v~ n t~ 3 , where v = r/t and n is 
a constant > 5. If the ejecta interact with a shallow power law medium p oc r~ s with s < 3, there 
exist self-similar solutions for the shocked flow (Chevalier 1982a). The problem contains only two 
independent dimensional parameters: the coefficients of the supernova density and the external 
density. In the exponential case there is an additional parameter, the velocity scale height v e , so 
that the problem is not self-similar. 

The exponential profile is an oversimplification of the density profile, as can be seen in fig. 
1 of DC98. In a thermonuclear explosion, the elements of gas receive most of their energy from 
the burning of that gas, as opposed to a shock wave as in the case of core collapse supernovae. 
The expansion of these heated elements is what leads to the approximately exponential density 
distribution. If there is incomplete burning of the material, the internal energy per unit mass of the 
gas is lower that than in neighboring regions where complete burning has taken place. The gas is 
expected to evolve toward pressure equilibrium, compressing the incompletely burned region. The 
magnitude of the effect can be estimated by considering the relative energy densities produced by 
different burning processes, and allowing the gas with the greater burning to adiabatically expand 
to come into pressure equilibrium with the other gas. The difference in burning a C/O mixture to 
Si or to Ni does not lead to a large density difference: < 40% higher density in the Si gas. We will 
find that this is not sufficient to explain the clumping indicated in Tycho's remnant. If some C or 
O is left unburned, a greater degree of compression is expected, but this typically applies to only 
the outermost layers of the supernova. This effect can apparently be seen in models (e.g., the W7 
model density distribution shown in Fig. 1 of Branch et al. 1985). 

On a longer timescale, the density structure can be changed by the Ni bubble effect. The 
occurrence of SN 1987A brought the realization that the energy release by radioactive 56 Ni and 56 Co 
after the initial explosion could affect the density and composition structure within the supernova 
(e.g., Woosley 1988). Li, McCray, k Sunyaev (1993) analyzed the Fe/Ni/Co lines from SN 1987A 
and showed that the filling factor of Fe with velocity < 2500 km s _1 was > 0.3 although the 
mass fraction of this element was only ~ 0.01 in the same region. They attributed this effect to 
the expansion of gas heated by radioactive power deposition. Basko (1994) made more detailed 
calculations of the expansion effect. Li et al. (1993) noted that the structure should show up in 
young supernova remnants and briefly mentioned the Crab Nebula. Chevalier (2000) contrasted 
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the operation of the effect in Type la and in core collapse supernovae, and we follow his discussion 
of the Type la events here. 

We examine the expansion of a spherical volume of radioactive gas, following Basko (1994). To 
get a crude estimate of the effect, we assume that immediately after the explosion, the supernova 
gas is freely expanding and has constant density throughout. The velocity at the edge of a spherical 
clump of 56 Ni with mass Mni is 

U = 7.7 x 103 (J**-) V3 (-*-) V2 (4-) ~ 5/6 km S -K (5) 
VO.5M ; Vl0 51 er gs ; \M ch ) y > 

The deposition of radioactive energy in the 56 Ni gas leads to expansion so that the velocity of the 
outer edge of the region becomes Uoo . Conservation of energy leads to the equation 

X 5 -X 2 = §, (6) 

where X = Uqo/Uq and Q is the energy per gram deposited by radioactivity. The assumption of full 
deposition of energy in the gas implies Q = 3.69 x 10 16 ergs g" 1 for 56 Ni^ 56 Co and Q = 7.87 x 10 16 
ergs g _1 for 56 Co— > 56 Fe. The expansion is limited by the diffusion of either the 7-rays or the 
photons out of the high pressure bubble so that the full values of Q listed above may not be 
achieved. The bubble expansion is expected to drive a shell in the incompletely burned matter on 
the outside. Basko (1994) showed that the shell is subject to the Rayleigh- Taylor instability, but 
that the instability does not have much effect on the shell expansion. The shell may break into 
clumps. Equation (6) shows that if Uq is large, the expansion of the Ni bubble is relatively small 
and the factor increase in the filling factor (A 3 ) is small. Clumps can still be created in the swept 
up gas. On the other hand, small regions of 56 Ni have the potential to substantially increase their 
filling factor in the supernova gas. 

The reference values used in equation (5) are those thought to be typical of Type la supernovae. 
These supernovae typically reach maximum light on a timescale of 15 days. This is the timescale on 
which diffusion of the internal radiative energy becomes effective, which has the effect of smoothing 
pressure gradients and stopping Ni bubble expansion. The mean life of 56 Ni is 8.8 days and that 
of 56 Co is 114 days, so we estimate that the 56 Ni radioactive energy is deposited in the radioactive 
gas, but not that from 56 Co. After the expansion due to the radioactivity, the outer velocity of the 
Fe-rich region is 8, 400 km s _1 and about 0.13 M of matter has been swept into a shell around the 
Fe-rich region. These numbers are approximate and will change for the density profile of a specific 
supernova. 

As discussed in § 1, Tycho's remnant shows evidence for two knots that have moved with 
undecelerated motion. At a distance of 2.5 kpc, the velocities of the knots are about 8, 300 km s _1 
(Hughes 1997). This velocity is consistent with that expected in a clumpy shell formed by the 
expansion of radioactive gas. However, Fe is expected in a clump only if it formed as non-radioactive 
Fe. The Fe clump may have formed as 54 Fe. X-ray spectra of Tycho's remnant have been obtained 
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with the result that, in general, Fe line emitting gas is at a higher temperature and lower density 
than Si line emitting gas (Hwang, Hughes, & Petre 1998 and references therein). This is also 
expected as a result of the Ni bubble effect. 

These considerations show that there is some evidence in Type la supernova remnants for 
clumping surrounding a large Ni bubble. The agreement between the knot velocities in Tycho's 
remnant and those expected for the Ni bubble shell is probably fortuitous. There is some evidence 
for both Tycho and SN 1006 that they were subluminous supernovae (e.g., van den Bergh 1993; 
Schaefer 1996), which would imply a smaller amount of 56 Ni. For the subluminous Type la SN 
1991bg, Mazzali et al. (1997, and references therein) estimate the synthesis of 0.07 M of 56 Ni 
extending to a velocity of 5,000 km s" 1 . This mass is consistent with the upper limit on the Fe 
mass (< 0.16 M ) found by Hamilton et al. (1997) for SN 1006. 

The density contrast that is created by the nickel bubble is not well determined. In the 
calculations of Basko (1994), the density contrast at the edge of the bubble is a factor ~ 5 above 
the surroundings. However, in the case of a Type la supernova with a relatively much greater 
amount of radioactive Ni, clumps of non-radioactive material may be immersed in the radioactive 
Ni and subjected to a greater compression (Hoflich 2000). We have considered a simple model for 
this process in order to estimate the compression. We assume that the expanding supernova gas 
has a constant density and that the free expansion during the first day is negligibly affected by the 
addition of energy, as above. The energy equation can be written as 

1 dpV _ dQ dV 
7-1 dt ~ ~dt~ P ltt , ( ' 

where 7 is the adiabatic index, V is the volume, and dQ/dt is the heat addition from radioactivity. 
Over most of the evolution, the gas is radiation dominated and 7 = 4/3, although the thermal 
pressure may be significant in the initial phases. Provided t < 1 day, the power input from 56 Ni 
is approximately constant at 4.8 x 10 10 ergs g _1 s -1 . Consideration of eq. (7) shows that in 
non-radioactive gas the pressure falls as t~ A ; this also initially applies to the radioactive gas, but 
it switches to t~ 2 evolution (at t ~ 2 x 10 3 s for typical parameters). If the radioactive gas can 
keep non-radioactive gas at its pressure, the compression of the non-radioactive gas increases as 
t 3 / 2 and can reach ~ 200 after a day. However, a typical clump cannot be maintained in pressure 
equilibrium at the high pressure and a compression wave is driven into the clump from the outside. 
When the radiative diffusion timescale becomes less than the age, there is the possibility that 
the gas is more compressible and larger density inhomogeneities are created. This effect has been 
found in calculations of shock breakout from a star in supernova models (e.g., Chevalier 1981b). 
An important aspect of the effect is that radiation pressure dominates the gas pressure, which is 
the case here. These issues will have to be examined in detailed computations. 

This discussion shows that it is possible that Type la supernovae do contain inhomogeneities 
which could be observed in a supernova remnant. The Ni bubble effect occurs on a timescale of 10 
days, or a scale of 10 15 cm for a velocity of 10 4 km s _1 . The structure resulting from the Ni bubble 
should become frozen into the ejecta at this early time. 
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3. METHOD 

We used the 2-D (two-dimensional) code ZEUS2D based on a finite difference scheme to carry 
out the hydrodynamical simulations (Stone & Norman 1992). The code uses an artificial viscosity 
to smooth shock transitions. We first ran 1-D simulations with an inner exponential density profile 
(ejecta) interacting with an outer constant density ambient medium, starting 2 weeks after the 
supernova explosion. The initial density distribution evolved into an intershock structure consisting 
of reverse-shocked ejecta in the inner region and a forward-shocked ambient medium on the outside, 
separated by a contact discontinuity (DC98). The 1-D intershock profile was then used to initiate 
2-D simulations. The inner ejecta gas was freely-expanding, using an inflow boundary condition. 
The outer gas was at rest. Gas pressure in the unshocked gas was unimportant. We neglected the 
effects of magnetism, heat conduction, and radiation. Magnetic fields may play a role in suppressing 
heat conduction and viscosity, but are not expected to be dynamically significant in the supernova 
remnant. The radiative cooling time for the optically thin gas at an age < 300 yrs is greater than 
the age of the remnant. We assumed that the energy losses through X-ray emission are small so 
that the dynamics of the flow are not affected by gas cooling. We used an adiabatic index 7 = 5/3. 

In 1-D runs, the starting dynamical age determines the initial ejecta-ISM interface position. 
Simulations with different initial ages showed that the evolution quickly converges within a few 
radial doubling times. The solution is unique, independent of when the simulation starts (see DC98). 
More generally, since the ejecta density profile only depends on its explosion mass and kinetic energy, 
the interaction with an ambient medium can be described by a set of scaling parameters R', V, T' 
using M, E, and p am (see Truelove & McKee 1999; DC98) 

/ w y (iLfV 3 PC (8) 

K^ParnJ \M ch ) 

V'=(^r) V2 « 8.45 x 10 3 ( T ^-] km s~\ (9) 



M ) ■ \M/M ch 



r=-^248E- r ^-j n--yr, (10) 

where no = /9 am /(2.34 x 10~ 24 gm cm~ 3 ) and p am is the ambient density. The dimensional vari- 
ables r, v, t can be expressed in terms of the nondimensional quantities: r' = r/R', t' = t/T', 
and v' = v/V. Nondimensional solutions can be conveniently returned to dimensional ones by 
re-scaling, and new dimensional solutions for different M, E, and p am can be calculated. For a par- 
ticular density distribution, one evolutionary sequence in the nondimensional variables represents 
all possible dimensional solutions. 

For the exponential density profile, the 1-D solution (Fig. 1) entered the Sedov- Taylor self- 
similar blast wave evolution with r oc i - 4 at t' S> 1. As the ejecta density decreased with time, the 
reverse shock, initially moving outward in the stellar frame, began to move inward at t' 2.5. It 
reached the stellar center at t' « 8. 
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4. INSTABILITIES FROM SMALL PERTURBATIONS 

We simulated the Rayleigh- Taylor instabilities in the intershock region in two-dimensional 
spherical polar coordinates assuming 0-symmetry. The 2-D numerical grid was initialized with the 
1-D solution at t' = 0.00054 (0.134 yrs for the standard parameters E$\ = 1, M = M c h, and uq = 1) 
across the angular domain, with perturbation seeds in density, pressure, and velocity placed between 
the contact discontinuity and the reverse shock wave. The grid was radially expanding, following 
the intershock boundaries until the reverse shock radius turned over to the stellar center. The 
evolution was tracked for five orders of magnitude in time. Simulations initiated earlier presented 
numerical problems because the higher density contrast across the contact discontinuity reduced 
the numerical time step determined by the Courant condition. The 2-D grid was moved based 
on the 1-D results. We used grid wiggling in the 1-D simulations, tracking the forward and the 
reverse shocks. In this process, the grid expands and contracts so that the shock locations are kept 
at the same grid numbers. We fit the evolution of the radii and velocities with a linear function 
and a logarithmic polynomial, respectively. We then applied the 1-D smooth fits on the 2-D grid 
boundaries, avoiding significant numerical noise coming from the wiggling scheme. The radial zone 
boundaries had a velocity distribution linearly varying with distance to the inner boundary, so that 
the zone spacing was kept uniform during the expansion. We used 500 zones in the radial direction 
and varied the number of angular zones; the radial zone spacing was 3 x 10~ 6 R' at t'=0.01, and 
increased to 1.6 x 10~ 3 R' at £'=1.73. 

We investigated the growth of instabilities with various grid resolutions and perturbations. 
To reduce the computational expense, only a fraction of a quadrant centered at 9 = 45° was 
used. The run with the finest resolution had 2000 angular zones in 1/2 of a quadrant. The radial 
boundary conditions were inflow and outflow at the inner and outer sides, respectively. The angular 
boundaries were reflecting. We generally used a spherical harmonic function for the perturbation, 

Y lm (9,4>) = Pr (cos 9)e m t (11) 

where m = and the harmonic I (perturbation mode) is even, considering the (^-symmetry and 
the reflection symmetry about the equator. The angular perturbation is essentially the associated 
Legendre polynomial P™, which has an increasing amplitude toward the polar axis. We applied 
perturbations such that the amplitude was between 1% and 50% near 9 = 45°. The perturbed 
width was less than 40% of the distance between the reverse shock and the contact discontinuity, 
about 2.5% of the whole intershock width. 

4.1. Evolution 

Fig. (8) illustrates the evolution of the instability perturbed by a 1%, I = 100 mode using 2000 
angular zones on 1/2 of a quadrant centered on 9 = 45°. Small perturbations in the density grow 
and form spikes that protrude from the contact discontinuity into the forward-shocked region. As 
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the remnant expands, the Rayleigh- Taylor instabilities continue to grow to form mushroom-shaped 
caps while the stems become narrower. Secondary fingers develop at the same time, so that there 
are more fingers than the initial perturbation mode suggests. After reaching their maximum extent, 
the fingers fall to their sides and interfere with neighboring fingers. The forward growth of the caps 
is blocked by the drag of the flow and the fingers do not extend farther from the dense shell at the 
contact discontinuity. 

In the later stages, the Kelvin-Helmholtz instability takes over at the mushroom caps, creating 
vortex rings. The relative motion of flows between a finger and its surroundings bends the stem 
and disrupts the flow. Since the effect of drag increases as the finger grows outward, the continual 
shedding of mass at the finger's top eventually leaves insuffient mass to overpower the drag of 
the countersteaming flow. The mushroom cap falls off to the side, and the remaining filaments 
are swept back. Vorticity develops in the less dense regions left by the original mushroom caps. 
Globally, the vortex rings gradually come to dominate over the spikes, and only the stems remain 
recognizable. 

The Rayleigh- Taylor instability consistently builds up fingers with long wavelengths (low /). 
The dynamical stages of evolution can thus be distinguished by the dominant mode. We found 
four stages of evolution. The first stage is characterized by linear growth of the Rayleigh- Taylor 
instability; the perturbations grow out from the initial seeds and evolve towards a mushroom shape. 
Second, nonlinear growth of secondary fingers among primary ones becomes important. The fingers 
multiply and become congested. Third, mutual interference among fingers increases the wavelength 
of the dominant mode; vortex rings develop as the Kelvin-Helmholtz instability becomes active, 
although new fingers still develop, as in a convective roll. The instability is fully developed at this 
point. In the final stage, the intershock density contrast is reduced as the expansion decelerates; 
vortex rings and all the other features in the flow are dispersed, and the instabilities fade. 

Figs. (3) and (4) show the angle-averaged density profiles with various resolutions at an early 
phase and at t' = 1.6 (approximately the present epoch for Tycho). Compared to the 1-D solution, 
the reverse shock front smears out somewhat because it is corrugated, and the density peak at the 
contact discontinuity is not present. The forward shock front is less affected. Various resolutions 
give basically the same profile, although finer grid resolutions reveal larger fluctuations. The late 
phase shows less density contrast and the decline of instabilities. 

We compare the exponential case to the power law case of n = 7 with s = (interstellar 
medium) and s = 2 (circumstellar wind) in Chevalier et al. (1992). In the power law case, the 
conditions are self-similar, so that the evolution is quasi-steady. Convection cells continually stir 
up a region near the contact discontinuity. The density structure remains qualitatively similar with 
time. The difference in the two n = 7 cases with different ambient media is that for the s = 2 case 
the density structure is characterized by groups of slender blobs with a pyramid outline, while for 
the s = case it is characterized by individual mushroom shapes. The reason is that a constant 
density surrounding medium causes a smaller density contrast across the the unstable region, and 
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a negative gradient of entropy. Thus, the s = 2 case allows the initial blobs to grow and stretch 
out into narrow fingers, contrary to the s = case in which the flow is blocked from propagating 
further and ends in a mushroom cap. 

The s = power law case studied by Chevalier et al. (1992) shows that a larger power law 
index n increases the dominant mode I and the evolution goes from quasi-steady convection to 
intermittent growth of slender fingers, approaching in a sense the general profile of the low-n, s = 2 
cases. This is due to the larger density contrast between the shocked ambient medium and the 
shocked ejecta at larger n, and to the fact that the gradients in the shocked ambient region become 
shallower. The exponential case exhibits evolutionary properties similar to the power law case with 
decreasing index with time, in that longer wavelengths are built up throughout the evolution after 
the initial growth stage. 

We show in Fig. (5) the approximate power law index n e j = — din p/dlnr in the unshocked 
ejecta just inside the reverse shock for the exponential case. In the 2-D simulations, the reverse shock 
surface is corrugated, but its position does not deviate significantly from its position in the 1-D case. 
The position of the contact discontinuity found by dimensional analysis is R c cx t( n_3 )/( ra_s ) for a 
power law ejecta density profile and an ambient medium p oc r~ s . The power law index for n = 5 
is the same as that for a Sedov- Taylor blast wave. For n < 5 the forward shock is still expected to 
expand with the Sedov law; however, the contact discontinuity expands more slowly than r oc t 0A 
and steadily distances itself from the forward shock in the self-similar frame (Chevalier 1982a). 
Self-similar solutions are no longer possible and the reverse shock starts to move inward towards 
the center. For the exponential case, we expect that the deceleration at the contact discontinuity 
is reduced when the reverse shock is in ejecta with d In p/dlnr < 5. The instabilities should then 
lose the deceleration that drives them. The index n starts to fall below 5 at t' > 0.325, while the 
turnover of the reverse shock actually occurs at t' ~ 2.5; the motion of the reverse shock takes time 
to respond to the changing density profile. 

Both the density contrast and acceleration decline rapidly with time in the exponential case; 
the decline of the instabilities is an inevitable consequence. At the end of our simulation (t 1 = 2.16), 
there is little continuing development of the instability. 

4.2. Variation of the Perturbation and the Resolution 

The evolution of the instability is insensitive to the initial linear perturbation. We have applied 
the initial perturbation in various areas: between the contact discontinuity and the forward shock, 
ahead the forward shock, and inside the reverse shock. After the growth stage, the sizes of the 
unstable regions are similar in all the cases. Perturbations placed at the shocked ejecta near the 
contact discontinuity are the most effective in exciting the instability. 

On a low resolution grid of 500 x 200 zones on 1/2 of a quadrant, a 10% perturbation with 
/ = 100 initially has 12 Rayleigh- Taylor fingers, which subsequently double in number. The I = 50 
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mode generates three secondary fingers between every two primary fingers, resulting in 24 fingers 
like the I = 100 mode (Fig. 8). There appears to be a preferred mode shortly after the linear 
growth, independent of the initial conditions. However, with finer resolution more details between 
the primary fingers appear and the presence of a dominant mode becomes vague. More highly 
resolved fingers appear to be more slender and sharp, and the overall flow patterns thus appear 
different. On finer grids, the flow pattern approaches the pyramid shape of an n = 7, s = 2 power 
law density case; mushroom structures become congested and twisted at the base near the contact 
discontinuity. Toward the end of the evolution, a preferred mode in increasing wavelength is built 
up independent of grid coarseness or the initial conditions. The presence of the late-stage dominant 
mode can be understood as the exponential model having a continually decreasing power law index, 
although limited resolution can also fabricate a longer wavelength dominant mode. Our results are 
consistent with the studies of Rayleigh- Taylor instability growth in a stellar explosion by Fryxell, 
Miiller and Arnett (1991), who used a PPM (piecewise parabolic method) code and cylindrical 
coordinates. 

The morphology of the fingers on our higher resolution runs is similar to that of Kane, Drake, 
& Remington (1999) using the PPM code PROMETHEUS. The results of Chevalier et al. (1992) 
using the piecewise parabolic code VH-1, Jun & Norman (1996a) using ZEUS2D, and Fryxell et 
al. (1991) using PROMETHEUS appear to have lower resolution than our results and do not show 
the narrow structure that we find. The high resolution calculation of Miiller, Fryxell, & Arnett 
(1991) is closer to the morphology that we obtain. Higher resolution calculations introduce shorter 
wavelength numerical noise and allow small perturbations to develop. Small scale structures have 
a faster growth rate of the Rayleigh- Taylor instability in the linear regime, and this can be seen in 
our simulations. Fig. (8) shows that the initial growth of Rayleigh- Taylor fingers leads to fingers 
with a larger extent at high resolution, presumably because the smaller cross section gives less 
drag. However, this trend may begin to reverse at the highest resolution because the stems of the 
Rayleigh- Taylor fingers become more unstable (see t' = 0.003 in Fig. 8). As mentioned above, 
the extent of the unstable region shows little dependence on resolution in the fully developed 
regime. On finer grids, the fingers in the nonlinear regime show an increasing tendency to bend, 
which prompts stronger interaction among fingers and produces more complex structures. As the 
resolution increases, the mushroom caps break up into filaments; this trend continues at our highest 
resolution. The clear mushroom shapes found in lower resolution calculations may be an artifact 
of the limited resolution. 

In independent work that was concurrent with our study, Dwarkadas (2000) simulated the 
instabilities resulting from the interaction of supernova ejecta (with an exponential profile) with 
a surrounding medium. He presented the evolution of the instability for a case in which it grew 
from numerical noise. This perturbation is much smaller than the cases that we considered and 
the instability remained in the linear regime until the late damping stage of the instability. When 
he initiated the simulations with a 2 % perturbation, his results are consistent with ours. 

We have also simulated the evolution for an initial perturbation with 100% amplitude. After 
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a transition phase during which the intershock region became more than twice as broad as in the 
standard case, the evolution returned to the standard case. As found by Kane, Drake, & Remington 
(1999) in simulations with a power law supernova density profile, the standard evolutionary track 
is robust. 



5. EVOLUTION OF CLUMPS 

As discussed in § 1, recent X-ray (Hughes 1997; Hwang & Gotthelf 1997) and radio (Velazquez 
et al. 1998) observations of Tycho's remnant show two knots juxtaposed near the eastern edge. 
The outer outline of radio (Velazquez et al. 1998) and X-ray (Vancura et al. 1995) emission, 
which presumably defines the forward shock front, shows an outward protrusion surrounding the 
knots. The properties of undecelerated motion (Hughes 1997) and an outer shock protrusion 
cannot be explained by instabilities generated by linear pertubations. The Rayleigh- Taylor fingers 
from linear perturbations are substantially decelerated and do not affect the outer shock front. 
Although our computations are two-dimensional, we believe that going to three dimensions would 
not change this conclusion. Rayleigh- Taylor fingers are known to show somewhat greater growth 
in three dimensions during the initial growth phase (e.g., Kane et al. 2000 and references therein). 
However, in computations with power law supernova ejecta, Jun & Norman (1996b) find that the 
thickness of the mixing layer is very similar in two and three dimensions when the instability is 
fully developed (see their figs. 2 and 3). The instability is probably limited by the properties of 
the interaction region. Jun & Norman (1996b) do find that the stellar ejecta show smaller scale 
structure in three dimensions. 

We believe that the origin of the knots is best explained by nonlinear clumps of supernova ejecta 
expanding into the intershock region. We thus simulated the hydrodynamic evolution of ejecta 
clumps interacting with the intershock gas. We represented clumps as denser spheres superposed 
on the smooth exponential density profile. In the simulations, we initialized a clump at the polar 
angle 9 = 45° on 2-D computations like those described in the previous section. In two dimensions, 
the sphere is in fact a 3-D torus around the polar axis. A 3-D sphere can be simulated by placing 
the clump on the polar axis 9 = on a 2-D grid; however, there are then singularity problems on 
the symmetry axis of the grid, which lead to long Rayleigh- Taylor fingers on the axis due to the 
instabilities discussed in the previous section. Comparisons of the hydrodynamics of ISM clouds in 
spherical and toroidal morphologies indicate that geometry contributes little to the interaction (Jun 
& Jones 1999). To examine the effect of geometry, we carried out clump interaction simulations 
with clumps at 9 = 15° and 9 = 80°; the clump evolution did not substantially change. 

Under the same boundary conditions as used for the instability computations, the clump was 
included in the freely expanding ejecta. The clump first ran into the reverse shock and then moved 
forward through the intershock region. The initial clump as well as the ejecta were cold compared 
to the high temperature shocked ejecta. The gas pressure was unimportant, so that the size and 
the density contrast of the clump remained unchanged until the clump expanded into the reverse 
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shock front. 

The basic physics of the interaction is similar to that for the interaction of nonradiative blast 
waves with an interstellar cloud (Klein et al. 1994). As the shock wave moves past the cloud, a 
reflected wave moves back into the shocked intercloud gas and creates a bow shock. The incident 
shock also creates a transmitted 'cloud shock' that moves into the cloud and crushes it. The cloud 
shock propagates with a velocity v c « (pi/ 'p c ) l l 2 v s , where pi, p c , and v s are the density of the 
intercloud medium, the density of the cloud, and the shock velocity in the intercloud medium, re- 
spectively. A larger density contrast between the cloud and the intercloud medium causes a smaller 
velocity ratio of v c /v s , which helps the development of a shear flow, and consequently the Kelvin- 
Helmholtz instability at the cloud-intercloud interface. The Richtmyer-Meshkov instability powered 
by the impulsive acceleration takes place due to the impact of the shock wave on the upstream side 
of the cloud-intercloud surface. When the cloud shock exits the cloud, a rarefaction wave moves 
back into the cloud and causes expansion. This acceleration leads to the Rayleigh- Taylor instability 
on the upstream side of the cloud. The combined instabilities lead to the destruction of the cloud 
on a timescale of several times of the cloud-crushing time, t cc = R c /v c . 

The clump-young supernova remnant interaction is more complicated than the cloud-blast 
wave interaction because of the structure of the intershock structure, including the pre-existing 
large scale instability. The most important factor measuring the crushing strength of a cloud-shock 
interaction is the impact area per unit mass, which is inversely proportional to the initial density 
contrast, % = p c /pi- We have explored the evolution of a single clump, varying its size, density 
contrast, and initial impact time with the reverse shock. We examined four initial impact times, 
with the clump density contrast ranging between 3 and 100, and the clump size ao (as a fraction 
of the intershock width) below 1/3 (Table 1). We estimate the size of the Tycho X-ray knots to be 
about 1/5 the intershock width, although this is uncertain because the clumps are poorly resolved 
and we see only shocked parts of the clumps. For a size below 30% of the intershock width at 
the corresponding interaction time, we found that (a) an early clump initiated at t' = 0.011 must 
have a density contrast x > 10 to cause a protrusion on the remnant outline; a clump with X = 30 
causes a protrusion but the effect would have subsided long before the present epoch for Tycho; 
(b) a clump initiated at t' = 0.22 with % = 50 has a bulging effect on the remnant outline, but 
the outline returns to spherical symmetry by t' = 1; (c) a later clump initiated at t' = 0.86 with 
X = 100 is not dense enough to reach the forward shock at present epoch t' = 1.6 for Tycho. A light 
clump with x = 3 develops vorticity when encountering the reverse shock and is quickly destroyed. 

Figs. (8) and (8) show the interaction with a single clump for x = 100 and two different 
interaction times. After passing through the reverse shock, the cloud gradually became flattened 
and curved like a crescent. Material streamed out from the horns of the crescent; the ram pressure 
difference between the axis and the side of the cloud drove the mass loss. The pressure near the 
cloud axis was higher because of the additional ram pressure on the front face of the clump. At 
the rear of the clump the flow became turbulent and left a trail of vorticies. The clump expanded 
laterally as it approached the forward shock. The front of the clump snowplowed the material 
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ahead. The forward shock wave was distorted and a bulge formed. The material on the shock front 
became distributed into two lumps that gradually receded from each other. As the expansion of the 
supernova remnant continued, the bulge as well as the crescent clump immediately behind it lost 
their identities. Even a clump with x = 100 could not penetrate the forward shock when initiated 
at a later time (Fig. 8), but snowplowed the material ahead and caused a protrusion. 

The velocity of a shocked clump is determined by the drag of the surrounding material. A 
shocked clump becomes comoving with the postshock flow in approximately the time that the cloud 
sweeps up a column density of intercloud gas equal to its initial column density. The drag time is 
proportional to t^ rag ~ X^^co so that a denser clump travels faster in the remnant. The lateral 
expansion of a shocked clump into a cresent shape increases the cross section area and the drag, 
significantly decelerating the clump. 

For an exponential density profile of the ejecta, the clump evolution depends on the time of 
clump interaction because of the evolution of the supernova remnant. An earlier clump requires less 
density contrast to have a protrusion effect on the remnant outline (Fig. 8) . Early interactions give 
the clump a crescent shape without the action of the Rayleigh- Taylor instability on the clump's front 
because the higher flow velocity in the remnant delays the rarefaction wave traveling back upstream 
as the cloud-shock exits the clump. Later clumps tend to develop instabilities on the clump surface. 
The exponential model gives a larger velocity difference and a larger density difference between the 
ejecta and the reverse shock front at earlier epochs. The early clumps rapidly move into a lower 
density medium and are thus more robust, for a given value of x- 

For n = 0.8 cm" 3 , E = 10 51 ergs, and M = 1.4 M (DC98), Tycho's remnant has t' = 1.6 
at its present age of 427 years. If this model applies and the knots are in undecelerated motion, 
they must have crossed the reverse shock front at t' « 1.2. This places Tycho in the regime where 
X > 100 is required. 

In addition to the isolated knots, we considered groups of clumps. The Ni bubble process is 
likely to create a shell of clumpy ejecta and, as discussed in § 6, the observations of ejecta in Tycho's 
remnant suggest that widespread dumpiness may be present. In Fig. (8) we show two simulations 
with multiple clumps. In the top one, there is one band of clumps and in the second, there are two 
bands. The evolution of each forward clump is similar to that in the isolated clump case, but they 
combine to push the forward shock wave to larger radius. The ejecta move out to a relatively large 
radius, as may be required in Tycho's remnant. 

6. THE REMNANT OF TYCHO'S SUPERNOVA 

As discussed in § 1, Tycho's supernova was likely to be of Type la. The expansion rate of 
the remnant has been measured in the radio, giving an expansion parameter m = 0.47 ± 0.05 
(Strom, Shaver, & Goss 1982), 0.462 ±0.024 (Tan k Gull 1985), and 0.471 ±0.028 (Reynoso et al. 
1997). The optical filament yields an expansion parameter of 0.39 ±0.01 (Kamper & van den Bergh 
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1978), close to the adiabatic blast wave case. The slight discrepancy in these measurements can be 
explained by sampling variance; the optical measurement samples only the densest regions where 
the emissivity is the highest. HI observations do indicate a higher surrounding density to the NE 
(Reynoso et al. 1999). These radio expansion values indicate a global pre-Sedov stage and exclude 
a surrounding circumstellar medium, which would give an expansion parameter ~ 0.7 (Chevalier 
1982b; DC98). The exponential model is successful in reproducing the observed shock position and 
deceleration of Tycho's remnant. At Tycho's present age of 427 years (t' = 1.7) the model places 
the reverse shock and the forward shock at radii of 2.08 pc (r' = 0.98) and 3.09 pc (r' = 1.45), 
while they decelerate with expansion parameters of m = 0.15 and m = 0.47, respectively. 

The radio observations of Tycho by VLA by Reynoso et al. (1997) show that the NE and par- 
ticularly the adjacent SE parts protrude from the circular outline. The SE protrusion corresponds 
to the two X-ray knots. Near the sharp NE edge, Velazquez et al. (1998) note the presence of regu- 
larly spaced structure. The spatial regularity is interpreted by Velazquez et al. as Rayleigh- Taylor 
fingers still in the linear regime with a mode / ~ 30; they attribute the preferred wavenumber to 
the effects of viscosity. However, the importance of viscosity is uncertain because magnetic fields 
can inhibit its effect. We believe that a linear regime for the Rayleigh- Taylor instability is unlikely. 
Regardless of the supernova density model, the shock front has expanded by a sufficently large 
factor to allow the instabilities to evolve to their saturated phase if there are initial perturbations 
with amplitude > 1%. 

The regular structure observed in the radio image by Velaquez et al. (1998) is suggestive of 
a Rayleigh- Taylor instability, but there are problems with this interpretation in the context of our 
models. One is that at t' = 1.6, the region occupied by the unstable ejecta does not extend beyond 
85 % of the remnant radius. The structure observed by Velaquez et al. is typically at a larger 
radius and in fact extends to the edge of the remnant. Another problem is that the instability at 
this late stage does not have a clear regular finger structure in our models, but is dominated by 
vortex rings with an irregular distribution. The large radius of the structure suggests that it may 
be connected to perturbations in the surrounding medium as opposed to the ejecta. However, the 
outer outline of the shock front is smooth on the scale of the structure, as imaged at both optical 
and radio wavelengths. In addition, the X-ray image of Tycho's remnant (e.g., Vancura et al. 1995), 
which is dominated by line emission, shows strong emission at a large fractional radius. Because 
the image is likely to be dominated by emission from the ejecta (Vancura et al. 1995; Hwang & 
Gotthelf 1997), the evidence is that the outer structure is related to the ejecta. 

We believe that a solution to this problem is indicated by the presence of the X-ray knots 
with undecelerated motion. Two knots are clearly observed in one section, but it is likely that they 
represent a more widespread phenomenon. We suggest that there was a spherical shell of clumpy 
ejecta at a velocity ~ 6, 000 — 8, 000 km s _1 in the freely expanding ejecta that gives rise to the 
clumpy structure. Kane et al. (1999) previously suggested that the radio structure observed in 
Tycho near the forward shock wave may be related to nonlinear density variations in the freely 
expanding ejecta. There are reasons to believe that the clumps are restricted to a region of the 
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ejecta and not spread through all velocities. First, the deceleration parameter measured at X-ray 
wavelengths is significantly larger than that measured at radio or optical wavelengths (Hughes 
1997). If the clumps are widespread, they should come and go in a band that expands as does 
the emission at lower wavelengths. However, a band of clumps could give the observed difference 
in deceleration parameters. Another reason is that in SN 1006, the Fe II absorption line profiles 
are consistent with a distribution of gas at velocities < 5, 000 km s^ 1 like that in the W7 model 
(Fesen et al. 1988), which has an approximately exponential density distribution. The Si line shows 
absorption at redshifted velocities 3,000 — 7,000 km s" 1 , which can be attributed to moderately 
high velocity Si interacting with the surrounding medium. There is no blueshifted Si absorption. 
The data are consistent with the picture of smoothly distributed Fe on the inside, outside of which 
is Si gas that may be clumpy. 

The presence of a band of clumps is expected in the Ni bubble scenario for clump formation 
(see § 2). One prediction of this model is that the outer ejecta emission should be very clumpy when 
observed at high X-ray resolution (see Fig. (8). In addition to this component, there should be 
emission from a smaller radius, near the reverse shock front. This emission should be representative 
of lower density, higher temperature gas. We believe that the emission in the Fe K line, which is at 
a smaller radius (Hwang & Gotthelf 1997), is primarily from gas close to the reverse shock wave. 
This gas should have some structure because of the earlier instabilities, but much less than the 
outer clumpy ejecta. 

Our computations showed that in order for knots to survive to close to the outer shock wave, 
they must have an initial density contrast with respect to their surroundings of ~ 100. Hwang et 
al. (1998) in fact find that the ionization time (n e t, where t is the time since the gas was heated) is 
100 times larger for Si than for Fe, based on line emission from the entire remnant. The time since 
heating should be comparable for the two elements, so this result is consistent with the Si being 
primarily in shocked clumps (created outside the Ni bubble) and the Fe primarily in interclump gas 
(inside the Ni bubble). However, the situation is complicated by the fact that there is evidence for 
some Fe in knots, which we argue had a nonradiogenic origin (e.g., 54 Fe). In addition, we estimated 
in § 2 that ~ 0.13 M of Si could be swept into a shell around the Ni bubble, so there could be some 
Si emission from the ejecta that is not clumped. Investigation of these issues will require spatially 
resolved X-ray spectroscopy, as will be possible with the Chandra and XMM observatories. 



7. IMPLICATIONS FOR TYPE la SUPERNOVA MODELS 

One result from observations of the Type la remnants is a lower limit on the velocity of 
Si in the explosion. For Tycho's supernova, there is an undecelerated knot of Si moving at 
8, 300(D/2.5 kpc) km s^ 1 , where D is the distance. In SN 1006, there is Si freely expanding 
with a minimum velocity in the range 5,600 — 7,000 km s _1 . In both of these cases, the veloci- 
ties appear to be less than the minimum velocity of the Si-dominant layer in the W7 deflagration 
model, about 10,000 km s _1 (Thielemann et al. 1986). However, if these are the remnants of 
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subluminous supernovae, the W7 model may not be the appropriate one. Also, the initial burning 
front is unlikely to be precisely spherically symmetric, so that some material that has not burned 
to 56 Ni may lag behind. 

Our most significant results involve the evidence for clumping in Type la supernovae. We have 
found that density contrasts > 100 are necessary to reproduce the large radius of X-ray emitting 
clumps and their velocities. A plausible mechanism for the formation of the clumps is the action of 
an expanding Ni bubble, although calculations of this process have not yet shown that such large 
compressions are feasible. We estimate that ~ 0.1 M of matter can be swept into clumps around 
the Ni bubble. For the two undecelerated X-ray emitting clumps in Tycho's remnant, Hughes 
(1997) estimates a mass of order 0.002 M Q for the Si + S clump and 0.0004 M for the Fe clump, 
so there can be many comparable clumps. The presence of Fe in a clump requires that it have a 
non-radioactive origin. The 1-dimensional calculations of nucleosynthesis associated with the W7 
Type la model (Thielemann et al. 1986; Iwamoto et al. 1999) show that 54 Fe is present along 
with Si and S in the region immediately outside of the 56 Ni region. However, the 54 Fe does not 
dominate the composition, unlike what is indicated by the observations of the Fe clump in Tycho. 
The central region of the W7 model is composed of stable 56 Fe and 54 Fe, but this gas is at too low 
a velocity to be compatible with the Tycho clump. Thus, the synthesis of fast material dominated 
by non-radioactive Fe remains to be demonstrated. 



8. CONCLUSIONS 

We used 2-D hydrodynamical simulations in spherical polar coordinates to investigate insta- 
bilities and dumpiness in a Type la supernova remnant. We adopted an exponential density profile 
for the smooth supernova ejecta and a constant density for the surrounding (interstellar) medium. 
The exponential density profile gives the best general approximation for the supernova ejecta in 
detailed explosion models, as opposed to a power law density profile for which self-similar solutions 
are available. The region between the forward shock in the interstellar medium and the reverse 
shock in the ejecta is Rayleigh- Taylor unstable. An initial perturbation grows to nonlinear strength, 
develops an increasing characteristic wavelength, and fades as the remnant evolves to the Sedov 
blast wave regime. The characteristic Rayleigh- Taylor mushroom caps are replaced by vortex rings 
from the Kelvin-Helmholtz instability, and turbulence appearing in vorticies remains dominant 
throughout the rest of the evolution. In these late phases, the depletion of the kinetic energy in 
the inner ejecta as the ejecta density profile continually weakens the acceleration needed to power 
the Rayleigh- Taylor instability. A sequence of simulations with increasing grid resolution shows 
the Rayleigh- Taylor fingers become increasingly slender and sharp and tend to bend, limited only 
by the grid resolution. 

The structure that is formed as a result of the Rayleigh- Taylor instability is decelerated and is 
confined to a region < 85% of the remnant radius for a remnant with the dynamical age of Tycho's 
remnant. These properties are not consistent with the knots observed in Tycho's remnant. We 
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thus considered a clump in the diffuse supernova ejecta with the aim of reproducing the properties 
of the the SE X-ray knots and the related protrusion on the forward shock outline. In the clump- 
shock interaction, the density contrast, clump size, and time of initiation of the interaction are the 
major factors. Ram pressure stripping gives rise to a core-plume structure and a strong Kelvin- 
Helmholtz instability in the downstream region. The Rayleigh- Taylor instability develops on the 
clump's upstream side, which facilitates fragmentation of the clump. The clump causes a bulge on 
the remnant outline as the ram pressure pushes material ahead. The bulge eventually disappears 
as the clump is fragmented and swept back. Punching through the forward shock does not occur 
for our similution with conditions like those in Tycho's remnant, even for the case of a clump with 
an initial density contrast of 100. Clump interaction in the early phases of the remnant evolution 
are more likely to distort the forward shock front because of the larger density contrast between 
the ejecta and the interstellar medium. 

The two X-ray clumps to the SE in Tycho are probably not alone. In fact, we find that the 
large radial position of most of the ejecta emission suggests that much of the emission is from 
clumps. If it were not in clumps, it could not extend to the observed radius. We thus suggest 
that the undecelerated knots are at the high column density end of a spectrum of knot properties. 
The required knot compression factor, up to ~ 100, is not a natural result of existing Type la 
supernova models. We suggest that the expansion of a radioactive Nickel bubble is responsible 
for the compression of the knots. The compression factor is larger than that found in existing 
calculations, but it is possible that radiative transfer effects during the final stages of the Ni bubble 
expansion can lead to enhanced compression. The Ni bubble interpretation has the implication 
that Fe clumps, which are observed, cannot have been synthesized as 56 Ni. They must have a 
non-radiogenic origin, such as 54 Fe. The synthesis of high velocity matter in which 54 Fe is the 
dominant component has yet to be shown in computations of nucleosynthesis. 

In addition to these suggested theoretical investigations, advances can be expected on the 
observational front. Spatially resolved X-ray spectra of Type la supernova remnants, as is possible 
with the Chandra and XMM observatories, will be valuable for showing whether our picture of 
ejecta clumps is valid. We predict that the outer ejecta emission in Tycho's remnant has a clumpier 
structure than the inner emission observed in the Fe K line. 
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Fig. 1. — One-dimensional evolution of the forward shock, the contact discontinuity, and the reverse 
shock radius with time. The farthest position to which the Rayleigh- Taylor fingers stretch in two- 
dimensional simulations is also plotted. The time and radius use normalized units given in the 
text. 
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Fig. 2. — evo.gif 

Series of density images illustrating the time evolution of the dynamical instability. The 
perturbation was imposed at t' = 0.00054 with 1% amplitude and an I = 100 mode of the spherical 
harmonic function. The grid has 500 radial zones by 2000 angular zones in 1/2 of a quadrant. 
The plots at early stages display the central 1000 angular zones in 1/4 of a quadrant, while the 
plots shown for later stages display 1/2 of a quadrant. The contour levels are exponentially spaced 
between the lowest value and the highest one sampled in the grid domain. 



-25- 



100 



10 



n i i i i i i i i r 



n i r 



V = 0.0027 




■,\ • 



■: \\ 

\ \\ 
\\\ 











500 90 
500 45 
1000 45 l 
2000 45 [ 



J I I I I I I I L_ 



J I I l_ 



0.0130 0.0135 0.0140 
Radius (r r ) 



0.0145 



Fig. 3. — Angle-averaged two-dimensional density distribution overplotted on the one-dimensional 
unperturbed solution in varying grid resolutions with I = 100 and 1% perturbing amplitude. The 
plot is at an early dynamical stage: t' = 0.0027. 
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Fig. 4. — Angle-averaged two-dimensional density distribution. The plot is at t' = 1.62, corre- 
sponding to the present age of Tycho's remnant. 
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Fig. 5. — Time evolution of the power law index of the ejecta density profile immediately below 
the reverse shock front from both the one-dimensional and two-dimensional runs. 
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Fig. 6. — resol.gif 

Series of density contours with varying resolutions for four epochs: t' = 0.0016, 0.0027, 0.011, 
and 0.11. The lowest resolution uses 200 angular zones on 1/2 of a quadrant and the highest uses 
2000 zones on 1/2 of a quadrant. The top two low resolution plots exhibit the same dominant 
modes on 1/2 of a quadrant, growing from the initial perturbation modes I = 100 and I = 50. The 
other plots are shown on 1/4 of a quadrant. All of the plots use a 10% perturbing amplitude. 

Fig. 7.— G0mplt.gif 

Snapshots of a clump in the ejecta expanding into the shocked region. The clump was initiated 
at t' = 0.207, with an initial density contrast x = 100 and a size ao = 1/8. The clump becomes 
crescent shaped and causes a protrusion in the forward shock front. The grid has 500 radial by 500 
angular zones on 1/2 of a quadrant. 

Fig. 8.— k2mplt.gif 

Snapshots of a clump with an initial density contrast x = 100 and a size ao = 1/6 in the ejecta 
expanding into the shocked region. The interaction was initiated at t' = 0.866. The Rayleigh- Taylor 
instability can be seen on the clump's downstream side. The grid has 500 radial by 500 angular 
zones on 1/2 of a quadrant. 

Fig. 9. — band.gif 

Snapshots of multiple clumps with a initial density contrast x = 50 and a size ao = 0.15 in the 
ejecta expanding into the shocked region. The interaction was initiated at t' = 0.866. The grid has 
500 radial by 500 angular zones on 1/2 of a quadrant. The top plot has one row of clumps and the 
bottom has two rows. 
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Table 1: Clump Simulation Parameters 



t' a 


o b 

He 


p c 
n rs 


Rfs d 


dR c 


0.011 


15.5 


0.042 


0.048 


0.006 


0.217 


2.03 


0.797 


0.958 


0.161 


0.866 


0.580 


1.70 


2.26 


0.478 


1.24 


0.315 


1.93 


2.70 


0.768 



"starting time of the clump-reverse shock interaction 

6 clump density for \ = 1, normalized to the unshocked ISM 

formalized radius of the reverse shock 

^normalized radius of the forward shock 

e normalized intershock width 
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